Cluster-size heterogeneity in the two-dimensional Ising model 
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We numerically investigate the heterogeneity in cluster sizes in the two-dimensional Ising model 
and verify its scaling form recently proposed in the context of percolation problems [Phys. Rev. 
E 84, 010101 (R) (2011)]. The scaling exponents obtained via the finite-size scaling analysis are 
shown to be consistent with theoretical values of the fractal dimension df and the Fisher exponent 
r for the cluster distribution. We also point out that strong finite-size effects exist due to the 
geometric nature of the cluster-size heterogeneity. 
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In studies of critical phenomena, the finite-size scaling (FSS) has proved to be an ex- 
tremely fruitful approach, owing to the recent development in computation It is through 
a crossover effect that this method connects numerical results to theoretical understanding 
so successfully. Although a typical simulation of a lattice system can only deal with a finite 
length scale L, one may regard L as a relevant scaling variable in the following way 
Under a scaling transformation by zooming factor b, i.e., L^ 1 — > bL^ 1 — > b 2 L~ x ■■•—>■ 6 n L _1 , 
the scale invariance of critical phenomena suggests that the singular part of the free energy 
/ will transform as 

/(tit, L" 1 ) = b- d f(b Vt u t , bL~ l ) = b- nd f(b nyt u t , b n L^), (1) 

where Ut is a thermal scaling variable and yt = \jv is its eigenvalue. The quantity v is 
called the correlation- length exponent since the correlation length £ behaves as £ ~ \A~ V 
with t = (T — T c )/T c , where T c denotes the critical temperature. The scaling variable u t is 
proportional to t/t near a fixed point of this renormalization-group (RG) transformation, 
where to is a nonuniversal scale specific to the system. This linear approximation remains 
valid until the system reaches a point where b nyt \t/to\ ~ 0(1). Comparing this with Eq. (JTj), 
we find that / is only weakly dependent on the first argument in this region and thus can 
be written as 

f{u ti L- 1 )r.\t/t Q \ dv f{\tlt Q \- v L- 1 ). (2) 

In other words, the crossover occurs when |t| ~ L^ 1 ^, and therefore, L usually enters an 
FSS form accompanied by the exponent —\jv. Recalling that £ ~ we see that the 

argument on the right-hand side of Eq. (j2J) means £/L, expressing the competition between 
£ and L. 

The cluster-size heterogeneity H suggested by Lee et al. has been devised in the context 
of recent debates on explosive percolation and defined as the number of distinct cluster 
sizes {3J]. For brevity, we will call it simply "the heterogeneity" henceforth, but it is to be 
noted that only sizes (or volumes) of clusters matter in calculation of H, irrespective of 
shapes of clusters. Noh et al. have shown that this quantity can be applied to ordinary 
percolation as well 4|: They have calculated H = H(p) with finite L's, where p means 
occupation probability, and found that H(p) has a peak at p = p*, which approaches the 
true percolation threshold p c as L increases. It is a typical signature to detect the percolation 
transition. One may well expect that the deviation of p* from the true p c , related to the 



thermal scaling variable in percolation, will be described as ~ L~ x l v for the same reason as 
explained above. However, Noh et al. [4j have revealed that the FSS form for H is obtained 
with another exponent ur than v. This exponent uh can be argued in the following way: 
Let us consider the probability distribution of cluster sizes s, which scales as 



P(s) ~ s~ T e 



s/s 



with the characteristic cluster size s c ~ Itl" 1 /"". Here t = p — p c in percolation problem and 
the scaling exponent r is called the Fisher exponent. If the system is off-critical, its cluster 
sizes will be simply found between 1 and s c so that H ~ s c ~ Itl" 1 /"". On the other hand, 
if \t\ <C 1, H is limited by the finite size L and its scaling relation becomes different: We 
first consider s n such that P(s n ) ~ s~ T ~ 1/N C , where N c is the total number of clusters. In 
fact, since the typical cluster size is sP(s) ~ 0(1), N c is comparable to the total number 
of points N. For a d- dimensional system, it is obvious that L is related to N by N ~ L d . 
Furthermore, since P{s n ) ~ 1/N C , there are few clusters above s n and their contribution to 
H will be much smaller than those with s < s n . It is likely to find at least one cluster for 
each s < s n since P(s < s n ) > l/N c , and thus it is plausible that H ~ s n ~ A^ 1 ^. To sum 
up, H is expected to have the following behavior 



H(t,L) -»■ 



L d l T if \t\ < 1 
kl- 1 / " otherwise. 



Note that it is not £ and L that actually compete here. The competition is rather between 
s c and s n , reflecting the limitation in observing large clusters due to the finite system size. 
The appropriate FSS form for describing H should be, therefore, 

H{t,L) = L d / T F H (\t\L 1 / u ») (3) 

with Uh = t /(da). This FSS form is readily supported by numerical results in percolation (4]. 

The ferromagnetic Ising model is one of the simplest and the most important models 
in statistical mechanics. When the dimensionality d is higher than unity, it undergoes an 
order- disorder transition at a critical temperature T c . Near the critical temperature T c , the 
magnetic order parameter m scales as a power-law form m ~ \tf and the corresponding FSS 
form is given as m = L~^/ u F m {\t\L l / v ) l|, |5j. Instead of this conventional order parameter, 



we will look into the system in terms of cluster statistics such as heterogeneity in the present 



study, where a cluster in the Ising model is defined as a set of nearest-neighboring spins with 
the same direction. We first recall that Suzuki jg| has conjectured that the fractal dimension 
of an Ising cluster is df = d—j3/v. For instance, for d = 2, this conjecture suggests df = 15/8 
by inserting (3 = 1/8 and v = 1 into the relation. However, this suggested fractal dimension 
is not consistent with numerical results [7|. The reason is that there are actually two 
contributions in forming a cluster: One is due to the correlation due to the spin interaction, 
while the other is a purely geometric contribution which survives even in the infinite-T 
limit js|. For example, one has a chance to find a giant cluster in the Ising model in the 
triangular lattice at T — > oo since each site has a spin state of either up or down randomly 
with probability 1/2, which coincides with the site-percolation threshold in the triangular 
lattice {9]. In order to separate the geometric effect from the correlation effect, one may 



ll|. These bonds are only 



consider bond-occupation probability pb inside a cluster 

ay 

to define connectivity between neighboring spins in the same direction and do not affect the 
spin interaction. Therefore, we are back to the original cluster statistics with ps = 1, which 
is of our primary concern. When ps = 0, on the other hand, the system reduces to the pure 
Ising model, with two relevant eigenvalues associated with t and h, respectively, where h is an 
external magnetic field. It means that we need think of the RG parameter space as (t,h,pB). 
Since pb is related to the random bond percolation, the p^-axis will have one unstable fixed 
point, separating two stable fixed points. A careful RG analysis shows that the stable fixed 
point at higher ps is actually a tricritical point ll|, ll2j. Suzuki's argument that df = 15/8 
is indeed true at the unstable fixed point describing Fortuin-Kasteleyn clusters [3], but one 
should note that the critical behavior is not given by this point, because it is the tricritical 
point that attracts the RG flow starting from pg = 1. Conformal invariance then predicts 
that the fractal dimension for geometrical clusters of the two-dimensional (2D) Ising model 



is df = 187/96 from tricritical 



numerical results in Refs. 
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16]. 



exponents 



12j. This prediction is well substantiated by 



Now we consider the cluster-size distribution in the square lattice at the critical temper- 
ature T c = 2/ln(l + \/2) ~ 2.27 in units of J/kg, where kg is the Boltzmann constant. 
This distribution follows a power law with exponent r in the thermodynamic limit, and the 
exponent r follows the relation [17| . 



r = 1 + d/df, 
4 



(4) 



which is derived from cluster statistics as in percolation. In other words, one can check the 
predictions of df by measuring r from the cluster-size distribution. 

We perform Monte Carlo (MC) simulations of the Ising model using the Metropolis and 
the Wolff algorithms [l| in 2D L x L square with under the periodic-boundary condition. 
For most simulations we use L = 20, 40, 60, 80, 100 and 160, but L = 320 is also used when 
required. We mostly use the Wolff algorithm for efficiency, and the Metropolis algorithm 
only for a consistency check. We start from a temperature much higher than T c , and then 
slowly decrease T, measuring equilibrium quantities at each temperature. All results are 
obtained from averages over 5 x 10 5 MC steps, after disregarding the first 5 x 10 5 MC 
steps for equilibration. We first choose a snapshot of a spin configuration in equilibrium 
and identify all the clusters in the system. Then we count how many different sizes of the 
clusters are found in spin-up and down directions, respectively. The cluster heterogeneity 
is then calculated as the sum of these two numbers, one in the spin-up direction, the other 
in the spin-down direction. Our brute-force approach in counting all the cluster sizes works 
as a main bottleneck in increasing the system size: For example, the CPU time spent for 
L = 160 roughly amounts to 4 x 10 2 hours. Although limited by such a practical difficulty, 
our results nevertheless nicely agree with the prediction in the conformal field theory (see 
Table I]). 



The cluster distribution in the 2D Ising model obeys the following form 1_7[ : 



n s = s- T F na (\t\s x ), (5) 

where n s is number of clusters with size s. This scaling form is valid near T c and with 
zero magnetic field. At T = T c , it leads to n s ~ s~ T . The scaling relation Eq. fll]) with 
df = 187/96 for the 2D Ising model gives us r = 379/187 ~ 2.03, which is shown to be 
consistent with our numerical results as shown in Fig. Ufa). For a more careful analysis, we 
allow the correction to the scaling form as n s ~ s~ r (l + As" A + . . .) with a constant A and 
a correction exponent A > and apply it for L = 160 and 320. We find that our numerical 
results are described sufficiently well by r = 2.01(2) as shown in Fig. QXb). However, more 
precise estimation of r is a formidable task, since it depends on the choice of the scaling 
region of s and also a further increase of L can alter the estimation made for smaller sizes. 
A better way is then to cross check with the outcomes of scaling relations as to be discussed 
below. 
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FIG. 1. (Color online) (a) Cluster size distributions n s for the Ising model in the L x L square 
lattices at the critical temperature T c . This figure is plotted in the form of \n(n s s T ) versus the 
cluster size s with the predicted exponent r = 1 + d/df = 379/187 ~ 2.03. The errorbars are of 
about the same sizes as the symbols for the scaling regime of s < 1000 and omitted for better 
visibility, (b) Estimation of the Fisher exponent r for L = 320 (see text for details). Log-binned n s 
are shown in the form of n s s 2 ' versus s (symbols) and the lines are results from the curve-fitting 
to n s s T ~ 1 + As~^ with r fixed to 2.03 and 1.99. It is clearly seen that the correct value of r is 
between 1.99 and and 2.03 and we thus conclude r = 2.01(2). 
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FIG. 2. (Color online) Characteristic cluster size s c in the Lx L square lattices at T c ~ 2.27( J/ks)- 
The data points are scaled by L~ 1 ^ ua ^ with varying l/(va) around 187/96 ~ 1.95, so that all the 
data points lie on a horizontal line within error bars. This estimate results in l/{va) = 1.95(2). 



The peaks at the tail part in Fig. [U(a) are due to giant clusters, exclusion of which 
yields the monotonically decreasing distribution without the peaks instead. This allows 
us to approximate F na (\t\s x ) as an exponentially decaying function e _S//Sc where s c is the 
characteristic cluster size similar to the one used in percolation {4]. We may also identify 
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FIG. 3. (Color online) (a) Heterogeneity H is drawn as a function of T for the 2D Ising model in the 
L x L square lattices. The dotted line indicates T c , and error bars are smaller than the data points. 

(b) We obtain Xu = 1.00(2) by scaling H(t = 0, L) as a function of L for the largest sizes at T c . 

(c) Scaling by using the conventional FSS ansatz in Eq. ([7]) with v = 1 and Xh = 374/379 ~ 0.99. 

(d) Scaling collapse using Eq. © with v H = 379/192 » 1.97 and X H = 374/379 » 0.99. 

s c with a peak position of n s near the tail in Fig. Ufa), and either way gives similar scaling 
behavior, s c ~ |t| _1 /°". In order to obtain the exponent a, we apply the standard technique 
of the FSS method to s c : 

s c = L^F Sc (\t\L^) . (6) 

The exponent a can be then obtained by analyzing simulation results with varying L at 
T = T c , by which we estimate the value of the exponent 1 / (ua) ps 1.95(2) (see Fig. [2]). Using 
the correlation-length exponent v = 1 of the 2D Ising model, 1/cr ~ 1.95(2) is numerically 
obtained. The scaling relation a — 1/ (vdf) allows us to check our numerical results with the 
predicted value 1/cr = 187/96 ~ 1.95 as shown in Fig.[2j Again the agreement is compelling, 
and our estimate favors this value over Suzuki's conjecture 15/8 = 1.875. 

We next carry out the FSS analysis of H in Fig. |3j First, we scale our numerical results 
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FIG. 4. Heterogeneity in the square lattices as a function of L, when each site is assigned a random 
spin state between up and down at T — > oo. Note that the horizontal axis is drawn on a logarithmic 
scale. Error bars are shown but smaller than the data points. 

[see Fig. E£a)] using a conventional FSS form: 



with a scaling exponent Xh- Although we can find the scaling exponent Xh ~ 1.00(2) 
by numerically observing H(t = 0, L) as a function of L [see Fig. EJ(b)], we cannot observe 
scaling collapse with v — 1 as shown in Fig. EJ^c). Consequently, the argument of Fh in the 
conventional FSS form Eq. (|7j) cannot be valid for H. This implies also that the correct 
scaling form should take into account the competition between cluster sizes as in percolation. 
By the same reasoning as above, a new scaling form is expected to be 



where u H = r/(da). Using r = 379/187 and 1/er = df = 187/96, we expect the exponent 
to be v H = 379/192 « 1.97, which indeed leads to a scaling collapse in a good quality 
[see Fig. E(d)]- In addition, our numerical estimate d/r = 1.00(2) contains the theoretical 
prediction 374/379 ~ 0.99 within the error bar. We have obtained basically the same 
scaling behavior for heterogeneity in the triangular lattice as well (not shown), confirming 
its universality in the 2D Ising model. Our results are summarized in Table U 

An interesting point in Fig. E^a) is that H keeps increasing as L grows even though T 
is far higher than T c . Recall that all the correlation due to spin interaction is destroyed in 
the infinite-T limit, where we are back to the percolation case. For the square lattice, the 

n 

site-percolation threshold is p^ uare = 0.59274602(4) [18[, while the infinite-T Ising model 



H(t,L) = L x »F H (tL^), 



(7) 



H{t,L) = L d/r F H (tL 1/uH ), 



(8) 
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TABLE I. Scaling exponents for heterogeneity of the 2D Ising model. The numerical values are 
obtained by Monte Carlo calculations in this work, whereas the analytic values are from tricritical 
exponents in the conformal field theory. 



df t d/r uh 

numerical 1.95(2) 2.01(2) 1.00(2) 1.96(5) 

analytic 187/96 379/187 374/379 379/192 

ps 1.95 ps 2.03 ps 0.99 ps 1.97 



corresponds to p = 1/2 due to the up-down symmetry. Since Noh et al. have argued 
that H ~ InL when \p — p c \ 3> L^ 1 ^" [J], we expect that H should have the logarithmic 
divergence for the infinite-T Ising model, which is confirmed by our numerical results (Fig.Hj). 
It is slower than the power-law divergence of L d l T at T c , so there will develop a peak at 
T = T c at a large L. The existence of such a peak should be true for the Ising model in 
the triangular lattice, too: As mentioned above, one can find a percolating phase in the 
infinite-T limit since p* rian g le = 1/2. The divergence at infinite T will be therefore described 
by 2D percolation exponents such as H ~ £i82/i87^ which j s slightly slower than the Ising 



case of H ~ L 374 / 379 at T = T c triangle = 4/ In 3 « 3.64 [19|. Returning back to the square 
lattice, since we have a good reason to believe the existence of a peak at T = T c for a large L, 
the monotonic shape of the scaling function in Fig. [3]^b) suggests that our observation might 
be still subject to corrections to scaling. We have indeed numerically found that the the 
logarithmic function H(L) is hardly distinguishable from linear increase when L < O(10 3 ), 
which explains the increase of H(L) in Fig. [3ja). 



le 2D Ising model 
4j, instead of the 



In summary, we have shown that the cluster-size heterogeneity H of t 
is scaled by the recently suggested FSS form using the exponent u H 
correlation-length exponent v. The finite-size effects in measuring this quantity are quite 
substantial, especially if compared to those in percolation. We have argued the main reason 
that H does not converge to a constant as the system size increases when T > T c but still has 
weak divergence due to its geometric nature. In spite of this difficulty, the scaling exponents 
at T = T c are in nice agreement with theoretical predictions (see Table [I]). This justifies the 
validity of this observable as well as the use of v H instead of v in the FSS analysis of the 
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Ising system. 
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